WAVE PROPAGATION IN WAVEGUIDES WITH RANDOM BOUNDARIES 
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Abstract. We give a detailed analysis of long range cumulative scattering effects from rough boundaries in waveguides. 
We assume small random fluctuations of the boundaries and obtain a quantitative statistical description of the wave field. The 
method of solution is based on coordinate changes that straighten the boundaries. The resulting problem is similar from the 
mathematical point of view to that of wave propagation in random waveguides with interior inhomogencities. We quantify the 
net effect of scattering at the random boundaries and show how it differs from that of scattering by internal inhomogencities. 
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^1 1. Introduction. Wc consider acoustic waves propagating in a waveguide with axis along the range 

Q ■ direction z. In general, the waveguide effect may be due to boundaries or the variation of the wave speed 

with cross-range, as described for example in [131 HH]- We consider here only the case of waves trapped by 
boundaries, and take for simplicity the case of two dimensional waveguides with cross-section T) given by a 
bounded interval of the cross-range x. The results extend to three dimensional waveguides with bounded, 
simply connected cross-section P C . 

The pressure field x, z) satisfies the wave equation 

< 



> 
oo 



p(t,x,z) = F{t,x,z) , (1.1) 



with wave speed c{x) and source excitation modeled by F(t,x,z). Since the equation is linear, it suffices to 
consider a point-like source located at {xo,z = 0) and emitting a pulse signal f{t), 

F{t,x,z)^ f{t)5{x-XQ)5{z). (1.2) 

Solutions for distributed sources are easily obtained by superposing the wave fields computed here. 

The boundaries of the waveguide are rough in the sense that they have small variations around the 
values .T = and x ^ X , on a length scale comparable to the wavelength. Explicitly, we let 

B{z)<x<T{z), where |B(z)| < X, |T(z) - X| < X, (1.3) 

and take either Dirichlet boundary conditions 

p{t,x,z) = {), for .T = B(z) and x^T{z), (1.4) 

X 

' or mixed, Dirichlet and Neumann conditions 



p{t,x = B{z),z)=Q, ■^pit,x = Tiz),z) = 0, (1.5) 

where n is the unit normal to the boundary x = T{z). 

The goal of the paper is to quantify the long range effect of scattering at the rough boundaries. More 
explicitly, to characterize in detail the statistics of the random field p{t,x, z). This is useful in sensor array 
imaging, for designing robust source or target localization methods, as shown recently in [3] in waveguides 
with internal inhomogencities. Examples of other applications are in long range secure communications and 
time reversal in shallow water or in tunnels [H [14] . 

The paper is organized as follows. We begin in section [2] with the case of ideal waveguides, with straight 
boundaries B{z) ~ and T{z) = X, where energy propagates via guided modes that do not interact with 
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each other. Rough, randomly perturbed boundaries are introduced in section [31 The wave speed is assumed 
to be known and dependent only on the cross-range. Randomly perturbed wave speeds due to internal 
inhomogcncitics are considered in detail in [T31 HH IH [HI [S] . Our approach in section [Hj uses changes of 
coordinates that straighten the randomly perturbed boundaries. We carry out the analysis in detail for the 
case of Dirichlet boundary conditions p. 41) in sections [Hj and [H and discuss the results in section [5] The 
extension to the mixed boundary conditions (jl.Sp is presented in section [5] We end in section [7| with a 
summary. 

Our approach based on changes of coordinates that straighten the boundary leads to a transformed prob- 
lem that is similar from the mathematical point of view to that in waveguides with interior inhomogencities, 
so we can use the techniques from [T31 [TH [U [51 [S] to obtain the long range statistical characterization of the 
wave field in section [31 However, the cumulative scattering effects of rough boundaries are different from 
those of internal inhomogencities, as described in section [5l We quantify these effects by estimating in a 
high frequency regime three important, mode dependent length scales: the scattering mean free path, which 
is the distance over which the modes lose coherence, the transport mean free path, which is the distance 
over which the waves forget the initial direction, and the equipartition distance, over which the energy is 
uniformly distributed among the modes, independently of the initial conditions at the source. We show 
that the random boundaries affect most strongly the high order modes, which lose coherence rapidly, that is 
they have a short scattering mean free path. Furthermore, these modes do not exchange efficiently energy 
with the other modes, so they have a longer transport mean free path. The lower order modes can travel 
much longer distances before they lose their coherence and remarkably, their scattering mean free path is 
similar to the transport mean free path and to the equipartition distance. That is to say, in waveguides 
with random boundaries, when the waves travel distances that exceed the scattering mean free path of the 
low order modes, not only all the modes are incoherent, but also the energy is uniformly distributed among 
them. At such distances the wave field has lost all information about the cross-range location of the source 
in the waveguide. These results can be contrasted with the situation with waveguides with interior random 
inhomogencities, in which the main mechanism for the loss of coherence of the fields is the exchange of energy 
between neighboring modes |13[ [T^ SI [HI [i] , so the scattering mean free paths and the transport mean free 
paths are similar for all the modes. The low order modes lose coherence much faster than in waveguides 
with random boundaries, and the equipartition distance is longer than the scattering mean free path of these 
modes. 

2. Ideal waveguides. Ideal waveguides have straight boundaries x = and x ~ X. Using separa- 
tion of variables, we write the wave field as a superposition of waveguide modes. A waveguide mode is a 
monochromatic wave P{t,x,z) = P{uj, x, z)e~^'^* with frequency w, where P(uj,x,z) satisfies the Hclmholtz 
equation 

[dl + dl + uj^/c^ix)]P{uj,x,z) ^0, 2 e M, a; e (0,X), (2.1) 

and either Dirichlet or mixed, Dirichlet and Neumann homogeneous boundary conditions. The operator 
+ ijP' l(?[x) with either of these conditions is self-adjoint in L^(0, AT), and its spectrum consists of an 
infinite number of discrete eigenvalues {Aj(a;)}j>i, assumed sorted in descending order. There is a finite 
number N(u)) of positive eigenvalues and an infinite number of negative eigenvalues. The eigenfunctions 
4>j{u),x) arc real and form an orthonormal set 

X 

dx<j>j{uj,x)(l)i{uj,x) = Sji , j,l>l, (2.2) 

where Sji is the Kronecker delta symbol. 

For example, in homogeneous waveguides with c{x) = Cq, and for the Dirichlet boundary conditions, the 
eigenfunctions and eigenvalues are 

(^) = /|sin (j^y A, (c.) = ( J) ' [{kX/n f - f] , j - 1, 2, . . . (2.3) 

and the number of propagating modes is N{uj) ~ \_kX /tt\ , where \_y\ is the integer part of y and k = lo/co 
is the homogeneous wavenumber. 
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To simplify the analysis, we assume that the source emits a pulse f{t) with Fourier transform 



J —oc 

supported in a frequency band in which the number of positive eigenvalues is fixed, so we can set N{uj) = N. 
We also assume that there is no zero eigenvalue, and that the eigenvalues are simple. The positive eigenvalues 
define the modal wavenumbers /3j(w) = ^/Xjjuj) of the forward and backward propagating modes 

P,{lu, X, z) = x)e±'ft■('^)^ J = 1, . . . , iV. 
The infinitely many remaining modes are evanescent 



with wavenumber Pjito) = y^—Xj{uj) . 

The wave field p(t, x, z) due to the source located at (a;o, 0) is given by the superposition of Pj{i-o, x, z), 



p{t,x,z)= I * 



l(0,oo)(^) 



J 2n £^7^ .^..^/M^ ' 



=N+1 VftM 



4-00 



The first term is supported at positive range, and it consists of forward going modes with amplitudes 
ttj^o/ \fWi E^nd evanescent modes with amplitudes ej,o/ \fWj- The second term is supported at negative range, 
and it consists of backward going and evanescent modes. The modes do not interact with each other and 
their amplitudes 



</'j(w,a;o), j = l,...,iV, 



are determined by the source excitation (11.21) . which gives the jump conditions at z = 0, 

p(cj, x, z = 0^) — p(ti;, .T, z = 0^) = , 

dzv(yi, x,z = 0+) - x,z ^ 0') = f{Lo)5{x - xq) . 



(2.4) 



(2.5) 



We show next how to use the solution in the ideal waveguides as a reference for defining the wave field 
in the case of randomly perturbed boundaries. 

3. Waveguides with randomly perturbed boundaries. We consider a randomly perturbed section 
of an ideal waveguide, over the range interval z € [0,i/e^]. There are no perturbations for z < and 
z > L/e^ . The domain of the perturbed section is denoted by 



where 



rf" = {(x, z) e M2, B{z) <x< T{z), 0< z < L/e^] , 



B{z) = eXii{z) , T{z) = X[l + sv{z)], e < 1. 



(3.1) 



(3.2) 



Here v and /i are independent, zero-mean, stationary and ergodic random processes in z, with covariance 
function 



'Ri,{z) ^¥.[v{z + s)v{s)] and 7e^(z) E[/i(z + s)^(s)]. 
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(3.3) 



We assume that and ^{z) are bounded, at least twice differentiable with bounded derivatives, and have 
enough decorrelatior0. The covariance functions are norniahzed so that TZ,^{0) and 7^.^(0) are of order one, 
and the magnitude of the fluctuations is scaled by the small, dimensionless parameter e. 

That the random fluctuations are confined to the range interval z G (0, i/e^), with L an order one length 
scale can be motivated as follows: By the hyperbolicity of the wave equation, we know that if we observe 
z) over a finite time window t € (0,T^), the wave field is affected only by the medium within a finite 
range from the source, directly proportional to the observation time T^. We wish to choose large 
enough, in order to capture the cumulative long range effects of scattering from the randomly perturbed 
boundaries. It turns out that these effects become significant over time scales of order so we take 

L"^ = i/e^. Furthermore, we are interested in the wave field to the right of the source, at positive range. 
We will see that the backscattered field is small and can be neglected when the conditions of the forward 
scattering approximation arc satisfied (see Subsection 14. 3p . Thus, the medium on the left of the source 
has negligible influence on p{t, x, z) for z > 0, and we may suppose that the boundaries are unperturbed at 
negative range. The analysis can be carried out when the conditions of the forward scattering approximation 
are not satisfied, at considerable complication of the calculations, as was done in |9j for waveguides with 
internal inhomogeneities. 

We assume here and in sections 2] and [S] the Dirichlet boundary conditions (|1.4p . The extensions to the 
mixed boundary conditions (jl.Sp are presented in section[6l The main result of this section is a closed system 
of random differential equations for the propagating waveguide modes, which describes the cumulative effect 
of scattering of the wave field by the random boundaries. We derive it in the following subsections and we 
analyze its solution in the long range limit in sectional 

3.1. Change of coordinates. Wc reformulate the problem in the randomly perturbed waveguide 
region by changing coordinates that straighten the boundaries, 

X = B{z) + [T{z) - B{z)] ^ , e e [0, X]. (3.4) 

We take this coordinate change because it is simple, but we show later, in section |4.4.2[ that the result is 
independent of the choice of the change of coordinates. In the new coordinate system, let 

u{t,^,z)=p[t,B{z) + [T{z)-B{z)]^,7}j , p{t,x,z) = u(t}^,^-^^^,z^ . (3.5) 

We obtain using the chain rule that the Fourier transform ^(a;,^, z) satisfies the equation 



1 + [(X - C)B' + £,T'] 



3 x'c,|s-*^ili^xc,.- 



{T -Bf ^ T -B 



2B'{T' - B') B" ^ 

{T- BY T -B ^ X 



T-B / T-B 



Xd(u 



+uj'^/^{B{z) + {T{z) - B{z))£_/X)u = , (3.6) 

for z e (0,i/e^) and ^ e {Q,X). Here the prime stands for the ^-derivative, and the boundary conditions at 
^ = and X are 

u{uj,0,z) ^u{uj,X,z) = 0. (3.7) 

Substituting definition p.2p of B{z) and T{z), and expanding the coefficients in p.6p in series of e, we obtain 

that 

(£0 + eCi + e^C2 + . . .) u[u, e, ^) = , (3.8) 



^Explicitly, they are 1/3-mixing processes, with 6 L^''^(R+), as stated in |15l 4.6.2]. 
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where 

C^=dl + dl+u^/c'{0 (3.9) 
is the unperturbed Helmholtz operator. The first and second order perturbation operators are given by 

£i + e£2 = q'{i, z)dl + M'{^, t z) , (3.10) 

with coefficient 

g^(e, z) = -2 [{X - Of^'iz) + ^i^'iz)] [1 - s {u{z) ii{z))] , (3.11) 

and differential operator 

M'{lo, i,z) = - |2 {v-^l)- 3e {v - fif - e [{X - Om' + ^i^f} 9| - 

{[{X Oa^" + e^"] [1 - e ('^ - A^)] - 2e (^^' - m') [(^Y - Oa^' + e^']} 5e + 



,2 



[(X - + en d^c--\0 + ^ [(^ - + ■ (3.12) 

The higher order terms are denoted by the dots in (|3.8p . and are negligible as e — >■ 0, over the long range 
scale L/e^ considered here. 



3.2. Wave decomposition and mode coupling. Equation p.Sp is not separable, and its solution is 
not a superposition of independent waveguide modes, as was the case in ideal waveguides. However, we have 
a perturbation problem, and we can use the completeness of the set of eigenfunctions (w, ^)}j>i in the 
ideal waveguide to decompose u in its propagating and evanescent components, 

N oo 

"(w,^,2:) = ^0j(w,^)%(cj,z) + ^ <f>j{uj,Ovj{uj,z). (3.13) 

3 = 1 j=N+l 

The propagating components Uj are decomposed further in the forward and backward going parts, with 
amplitudes aj{uj, z) and hj{u], z), 

u, = (S.e^ft- + 6,e-'ft^) , j = 1, . . . , iV. (3.14) 

This does not define uniquely the complex valued aj and hj , so we ask that they also satisfy 

5,% = i {a^e^P^^- - 6,e-*'3^-) , j = 1, . . . , TV. (3.15) 

This choice is motivated by the behavior of the solution in ideal waveguides, where the amplitudes are 
independent of range and completely determined by the source excitation. The expression (j3.13p of the wave 
field is similar to that in ideal waveguides, except that we have both forward and backward going modes, in 
addition to the evanescent modes, and the amplitudes of the modes are random functions of z. 

The modes are coupled due to scattering at the random boundaries, as described by the following system 
of random differential equations obtained by substituting p.l3p in p.8p , and using the orthogonality relation 
(|2.2p of the eigenfunctions, 

N . oo 

d.a, = [q,Sie^(^'-^^> + + ^e-'^^^ (gj, d.m + M^.di) + 0{s') , (3.16) 

1 = 1 Pj l=N+l 

N . oo 

d^b, ^ ^ \q, aze'(^'+^^> + C^6,e-(A-'^^)^j - ^e-'^^- (Q^=, d.di + M^, + 0{s^) . (3.17) 

1=1 Pj l=N+l 



The bar denotes complex conjugation, and the coefficients are defined below. The forward going amplitudes 
are determined at z = by the source excitation (recall (|2.4|) ') 



aj(w,0) = aj,o(w) , j^l,...,N, 



and we set 



0, j = l,...,N, 



(3.18) 



(3.19) 



because there is no incoming wave at the end of the domain. The equations for the amplitudes of the 
evanescent modes indexed hy j > N are 



N 



1=1 



-eJ2 {Qjid,di + M]idi) +Oie^), (3.20) 



l=N+l 



and we complement them with the decay condition at infinity 



lim Vj{uj,z) = 0, j > N. 

2— 5-±00 



(3.21) 



The coefficients 



z) = C]j^ {uj, z) + eCf (lu, z), for j > 1 and ; = 1, . . . , N, 



are defined by 



1 



X 



2y/pj{uj)l3i{uj) Jo 
1 



X 



2^/3, (w) A (w) 7o 
in terms of the linear differential operators 



and 



Ai = ~2{u - m)5| - 2jA [{X - Of^' + ^>^'] - [{X - Oa^" + e^"] + 

Bi = [z{v - m)' + [{X - + iA^] a| + 2iPi{v - ^) [{X - + e^^'] + 

{(j. - A.) [(X - OAi" + e^^"] + 2(i.' - //) [{X - + ^A} + 

'^[{X-Of^ + ^'^fdlc-^O- 



(3.22) 

(3.23) 
(3.24) 



(3.25) 



(3.26) 



We also let for j > 1 and / > N 



X 



1(2)/ 



X 



M^i{io,z) = I di<l,,{^,i)M'{uj,i,z)U^,0 ^ M'--l\uj,z)+eM];>{uj,z) 



(3.27) 



3.3. Analysis of the evanescent modes. We solve equations p.20[) with radiation conditions p.2ip 
in order to express tlie amplitude of the evanescent modes in terms of the amplitudes of the propagat- 
ing modes. The substitution of this expression in p.l6p - p.l7p gives a closed system of equations for the 
amplitudes of the propagating modes, as obtained in the next section. 

We begin by rewriting p.20p in short as 



l=N+l 



where 



and 



z) = gf\uj, z) + eg'/>{uj, z) + 0{e') , j > N, 



,(2) 



At) 



N 



r = 1,2 and j > N. 



1=1 



Using the Green's function Gj = e ^^1^1/(2/3^), satisfying 



d^Gj - /3|Gj = -d{z) , lim = , j> N. 



I 2:|— ^OO 



and integrating by parts, we get 



lil-s^)d]^iu;,z) 



2/3, (c.) 



(3.28) 



(3.29) 



(3.30) 



(3.31) 



(3.32) 



Here I is the identity and 5* is the linear integral operator 



-1 ^ />oo 

[^v],{u;,z) = -—- ^ / dse~^^^-^\^\{M^i-d.Q'^i){u;,z + s)M^,z + s) + 
- J2 / dse-^^^"^^'^sgn{s)Q^u,z + s)Mu^,z + s), 



(3.33) 



l=N+l ' 



acting on the infinite vector v ~ (utv+i, • ■ •) and returning an infinite vector with entries indexed by 

j, for j > N. The solvability of equation p.32[) follows from the following lemma proved in appendix El 

Lemma 3.1. Let he the space of square summable sequences o/L^(R) functions with linear weights, 
equipped with the norm 



\ 



j=N+l 



The linear operator "i! : Cn ^ defined component wise by p.33p is bounded. 
Thus, the inverse operator is 

(/-£*)-! =/ + £* + ..., 

and the solution of (|3.32l) is given by 



%(w,z) 



2/3, (w) 



ds e-^^(")l"lgj''(c^, z + s) + ©(e^) , 
7 



(3.34) 



Using definition p.30|) and the fact tliat tlie z derivatives of a; and hi are of order e, we get 



N 



N 



^Si(w,z)e-**M dse-^'^^^^'^-'^'^^^'C^^\uj,z + s) + Oie^). (3.35) 



We also need 



Wji(^,z) = dzVj{uj,z) 



(3.36) 



which we compute by taking a z derivative in p.28p and using the radiation condition Wj{ijj,z) — > as 
\z\ — >■ oo. The resulting equation is similar to p.32p 



h 2 



sgn(s)gj^(a;, z + s) + ^ Mj;(a;, z + s)vi{uj, z + s) 



l=N+l 



where we integrated by parts and introduced the linear integral operator 



[§w],{uj,z)^- / dse-^^(")l^lsgn(s)g|;(w,z + s)wi(w,z + s) 



(3.37) 



(3.38) 



l=N+l 



This operator is very similar to ^' and it is bounded, as follows from the proof in appendix [X] Moreover, 
substituting expression p.35p of vi in p.37p we obtain after a calculation that is similar to that in appendix 
|A]that the series in the index / is convergent. Therefore, the solution of p.37p is 



Wj{uj,z)^- I dse-ft(")l"lsgn(s)5^(w,z + s) + 0(£2) 



(3.39) 



and more explicitly, 

I ^ r°° 

d,VjiLo,z)^eJPjiuj)J2aii^,z)e'f^'" / ds e-^^^'^^\'^+"'''^^>sgnis)C^l\uj, z + s) + 

1 1 J — oo 



1=1 

N 



e^^j{uj)Ybi{uj,z)e-'^'' / dse-'^'('")l"l-"^'("'^sgn(s)cj/^(w,z + s) + 0(e2). (3.40) 

3.4. The closed system of equations for the propagating modes. The substitution of equations 
p.35p and p.40p in p.l6p and p.l7p gives the main result of this section: a closed system of differential 
equations for the propagating mode amplitudes. We write it in compact form using the 2N vector 



Xuj{z) 



a[uj, z) 
b{uj, z) 



(3.41) 



obtained by concatenating vectors a(aj, z) and b(a;, z) with components aj{u!, z) and 6j(w, z), for j — 1, . . . ^ N. 
We have 



with 2N X 2N complex matrices given in block form by 



H^(^) 



H^(z) HL°^(z) 

hL''(z) nt\z) 



(b), 



Gl^'jz) G^iz) 
gI'\z) Gi^\z) 



(3.42) 



(3.43) 



The entries of the blocks in are 



HL%iz) = ^C^\u.. z)e^^-^^^ , ifi^) (z) = , (3.44) 



and the entries of the blocks in G^j are 

r(l) 



°° IVv '((,1 7^ r°° 



"^31'^ / rfse-ft'H+'^A^A,sgn(s)4V(c.,z + s), (3.45) 



;'=A'+i 



°° '((,1 7^ r°° 

^^]}{u,z) 



l'=N+l 2 PjPl' J-oo 



,e-(A+ft)^y / dse-ft'l^l-'^^'^A,sgn(s)d,V(c.,z + s). (3.46) 

The coefficients in ^.44^ - ^.46^ are defined in terms of the random functions viz), fi{z), their derivatives, 
and the following integrals, 

1 '•^ 



- :rj^ / rfe</'j(0 h2a| + c.2e9^c-2(0] 0;(O , (3.47) 
^yPjPi Jo 



1 



c^^.iiio) = -j= \ di<i>,{i) {2dl^u?{X-i)d^c-\i)\ MO , (3.48) 
PjPi Jo 



1 



X 



d,M^)^--j== dUMOd^O, (3.49) 
^yPjPi JO 



d,,,i{^) = --^ I di{X~Om)dMO. (3.50) 



f-A' 

satisfying the symmetry relations 

c^j;(w) = Cpjj(w) , 
d^jiiuj) + d^^ij{uj) 



dM,,7(^)+dM,i,(w) = --7=%==. (3.51) 
2\/ft(^)/5i('^) 

We have from (|3?23| that 

= i^(z)c,ji(w) + K(z) + 2iA(w)zy'(z)]d,j,(L^) + 

A*(^)w(w) + + 2i/3i(w)/i'(z)]d^,j,(w) , (3.52) 

and from dS^Hl), ([XTTj) . ({31^ that 

— 2[v' {z)dyjv{uj) + ^' {z)d^ji'{LL>)] , 



M^]}{oj,z) 



iy{z)c„ji' (lo) + n{z)c^ji' (lo) + u" {z)d^^.ji' {lo) + ^j," (z)d^,^ji> (uj) . (3.53) 



4. The long range limit. In this section we use the system p.42p to quantify tlie cumulative scattering 
effects at the random boundaries. We begin with the long range scaling chosen so that these effects are 
significant. Then, wc explain why the backward going amplitudes arc small and can be neglected. This is 
the forward scattering approximation, which gives a closed system of random differential equations for the 

amplitudes {aj}j=i n- We use this system to derive the main result of the section, which says that the 

amplitudes {aj}j=i....^Ar converge in distribution as e — )■ to a diffusion Markov process, whose generator 
we compute explicitly. This allows us to calculate all the statistical moments of the wave field. 



4.1. Long range scaling. It is clear from p.4ip that since the right hand side is small, of order e, 
there is no net effect of scattering from the boundaries over ranges of order one. If we considered ranges of 
order 1/e, the resulting equation would have an order one right hand side given by H(^(z/e)Xtj(z/e), but 
this becomes ncghgible as well for e — >• 0, because the expectation of li^{z/e) is zero Chapter 6]. We 
need longer ranges, of order to see the effect of scattering from the randomly perturbed boundaries. 

Let then a^, bj be the rescaled amplitudes 

and obtain from p.42p that X^{z) ~ X^{z/e^) satisfies the equation 

dXUz) \^ ( ^e(^) ^ (^^^ xej^^ 0<z<L 



dz 



(4.1) 



(4.2) 



with boundary conditions 



(4.3) 



We can solve it using the complex valued, random propagator matrix P5(^;) G c2Afx2Af^ solution of the 
initial value problem 



'^P-(^) =iH^(^)p^(^) + G^(J)PS,(z) forz>0, andP^(0)=I. 



dz 



(4.4) 



The solution is 



XUz)=PUz) 



2o(w) 



and b'^{Lu,0) can be eliminated from the boundary identity 



ao(w) 



(4.5) 



Furthermore, it follows from the symmetry relations (|3.43p satisfied by the matrices and G;^ that the 
propagator has the block form 



P^"(^) 



P'Jiz) 



(4.6) 



where P5'°(z) and Plj^{z) are N x N complex matrices. The first block P^'" describes the coupling between 
different forward going modes, while P^''' describes the coupling between forward going and backward going 
modes. 

4.2. The diffusion approximation. The limit P^ as e — !• can be obtained and identified as a 
multi-dimensional diffusion process, meaning that the entries of the limit matrix satisfy a system of linear 
stochastic equations. This follows from the application of the diffusion approximation theorem proved in 
[18| . which applies to systems of the general form 



dX^iz 
dz 



-j^(x^z),y(^) ,^ 
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forz>0, and Af''(0) -Yo, (4.7) 



for a vector or matrix X'^{z) with real entries. The system is driven by a stationary, mean zero and mixing 
random process y{z). The functions T{x,y,T) and QiXiUi'T') are assumed at most linearly growing and 
smooth in and the dependence in r is periodic or almost periodic Section 6.5]. The function J-{x, y, t) 
must also be centered: For any fixed x f^i^d r. E[J^(x, 3^(0), r)] ~ 0. 

The diffusion approximation theorem states that as e — > 0, X'^{z) converges in distribution to the 
diffusion Markov process X(z) with generator C, acting on sufficiently smooth functions (p{x) a-s 

1 '•^ 



C^ix) = lim - / dr dzE [T{x. J^(0), r) • V^, [T{x, Vi^), r) ■ ^^Mx)]] + 

T-j-oo 1 Jq 







^ l\TE[gix,y{o),T)-vMx)] ■ (4.8) 



'0 

To apply it to the initial value problem (j4.4p for the complex 2N x 2N matrix P^(z), wc let X^{z) be 
the matrix obtained by concatenating the absolute values and phases of the entries in P^(2). The driving 
random process y is given by ijl{z),v{z) and their derivatives, which are stationary, mean zero and mixing 
by assumption. The expression of functions F and Q follows from (j4.4|) and the chain rule. The dependence 
on the fast variable r = z/e^ is in the arguments of cos and sin functions, the real and imaginary parts of 
the complex exponentials in (|3.44|) - (|3.46p . 

4.3. The forward scattering approximation. When we use the diffusion-approximation theorem 
in [18], we obtain that the limit entries of 'Pl;''{z) are coupled to the limit entries of 'Pl;°-{z) through the 
coefficients 

j'OC poo 

T^Af^j + A) - 2 / dzn,{z) cos[(/3j + ^i)z] , TZf^ifij + A) = 2 / dzTZf,{z) cos[(/3, + f3i)z] , 

JQ Jo 

for j, Z = 1, . . . , A''. Here TZi, and 7?.^ are the power spectral densities of the processes v and fi, the Fourier 
transform of their covariance functions. They are evaluated at the sum of the wavenumbers /3j + /3; because 

the phase factors present in the matrix Hi?''(z) are ±(/3j + /3i)z. The limit entries of Pfj°'{z) arc coupled to 
each other through the power spectral densities evaluated at the difference of the wavenumbers, TZ^{/3j — /?;) 
and TZ^iPj — Pi), for j,l = 1, . . . , A^, because the phase factors in the matrix H^"^ (z) are ±(/3j — Pi)z. Thus, if 
we assume that the power spectral densities are small at large frequencies, we may make the approximation 

^.(/3, +A)«0, n^iPj + /3i) ^ , for .7, / = !,..., TV, (4.9) 

which implies that we can neglect coupling between the forward and backward propagating modes as £ — > 0. 
The forward going modes remain coupled to each other, because at least some combinations of the indexes 
j, for instance those with \j — l\ = 1, give non-zero coupling coefficients TZi,{(3j — f3i) and ^^^{(Sj — Pi). 

Because the backward going mode amplitudes satisfy the homogeneous end condition 6^(cj, L) = 0, and 
because they are asymptotically uncoupled from {a^}j=i,...,Ar, we can set them to zero. This is the forward 
scattering approximation, where the forward propagating mode amplitudes satisfy the closed system 



da' 

dz e 



^H(f) (^) + G(,") (J) for z > 0, and a^co, z = 0) = %o(c^). (4.10) 



Remark 4.1. Note that the matrix hL"'' is not skew Hermitian, which implies that for a given e there 
is no conservation of energy of the forward propagating modes, over the randomly perturbed region, 

N N 

5]|a|(L)|VEl%°l'- 

This is due to the local exchange of energy between the propagating and evanescent modes. However, we will 
see that the energy of the forward propagating modes is conserved in the limit e ^ 0. 
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4.4. The coupled mode diffusion process. We now apply the diffusion approximation tlieorem to 
the system (j4.10p and obtain after a long calculation that we do not include for brevity, the main result of 
this section: 

Theorem 4.2. The complex mode amplitudes {a^(aj, z)}j=i^..._Ar converge in distribution as e — > to a 
diffusion Markov process process {aj(aj, z)}j=i^..._7v with generator C given below. 

Let us write the limit process as 

z) = P,{oj, J = 1, . . . ,iV, 

in terms of the power |aj P = Pj and the phase 9j. Then, we can express the infinitesimal generator C of the 
limit diffusion as the sum of two operators 



C = Cp + Ce 

The first is a partial differential operator in the powers 



N 



J, i = 1 



PiPj 



d d \ d 



dPj dPi J dPj 



( Pi - Pj 



dPi 



(4.11) 



(4.12) 



with matrix r^'^^(a;) of coefficients that are non-negative off the diagonal, and sum to zero in the rows 



(4.13) 



The off-diagonal entries are defined by the power spectral densities of the fluctuations and ^, and the 
derivatives of the eigenfunctions at the boundaries. 



4/3, (^) A (^) 



[ae0j(w,O)5e0z(^,O)]'^^[/3,(w) - Piiio)]} . 
The second partial differential operator is with respect to the phases 



(4.14) 



1 ^ 



Pi Pi . „ 



Pid9f^p,de]^^^^., 



de.dOi 



3,1=1 



(4.15) 



j.i = i 

3 it I 



with nonnegative coefficients 



(4.16) 



and 



x^ 



4ft- (w) A (w; 



[a^0,(c.,o)a^0,(^,o)]SM..K^)} , 



(4.17) 
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for j ^ I, where 



7i'j/(w) = 2 / dz sm[{l3j{uj) - I3i{uj))z]n^{z) , 



The diagonal part of r('''(a;) is defined by 



dz sin[(/3j(tj) - f3i{uj))z]Tlf_,{z) . 



(4.18) 
(4.19) 



(4.20) 

All the terms in the generator except for the last one in (|4.15p are due to the direct coupling of the propagating 

(4.21) 



modes. The coefficient Hj in the last term is 



with the first part due to the direct coupling of the propagating modes and given by 

N 



7^.(o) 



di 



4/3, 



(A + ft) Kji {f3? - ft') + 2d,,jic.,ji] 



^ terms, 



(4.22) 



with the abbreviation terms" for the similar contribution of the ^ process. The coupling via the evanescent 
modes determines the second term in (|4.2ip . and it is given by 



Q2\2 



i=N+i Wji3iW]+i3fy Jo 

oo 



ds e~P''nl{s) [{Pf - ft2) cos(ft s) - 2^j-A sin(;3js)] + 



dl.ijn'm + ^2T^7e.(o) 

ft Pi 



fi terms. 



(4.23) 



l=N+l 

4.4.1. Discussion. We now describe some properties of the diffusion process 2: 



1. Note that the coefficients of the partial derivatives in Pj of the infinitesimal generator C depend only 
on {Pi}i=i^...^N- This means that the mode powers {|a^(aj, z)P}j=i....,Af converge in distribution as 
e — >■ to the diffusion Markov process {|ft (a;, z)p = Pjito, z)}j=i,,,,^N, with generator £p. 

2. As we remarked before, the evanescent modes influence only the coefficient Kj(w) which appears 
in Ce but not in Cp. This means that the evanescent modes do not change the energy of the 
propagating modes in the limit e — )■ 0. They also do not affect the coupling of the modes of the 
limit process, because kj is in the diagonal part of (|4.15p . The only effect of the evanescent modes 
is a net dispersion (frequency dependent phase modulation) for each propagating mode. 

3. The generator C can also be written in the equivalent form [Sj Section 20.3] 



^ = I E rj-f (^) {A,^A,i + A,iA,i) + i ^ r<^\co)A,,A, 

3,1 = 1 j,l=l 

N 

E rjf -^«) + *E^^H^^-^-' 



(4.24) 



3,1 = 1 
3 1^1 



3 = 1 
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ill terms of the differential operators 



=^3^- = -AT, ■ (4.25) 



dai da 



Here the complex derivatives are defined in the standard way: ii z ~ x + iy, then dz = (l/2)(9j. —idy) 
anda^- (l/2)(a, + i9y). 

The coefficients of the second derivatives in (|4.24p are homogeneous of degree two. while the coeffi- 
cients of the first derivatives are homogeneous of degree one. This implies that we can write closed 
ordinary differential equations in the limit e — > for the moments of any order of {aJ}j=i,...,Ar. 

Because 

/ N \ 

" ~ =0, (4.26) 




we have conservation of energy of the limit diffusion process. More explicitly, the process is supported 
on the sphere in with center at zero and radius Ro determined by the initial condition 



N 



1=1 

Since C is not self-adjoint on the sphere, the process is not reversible. But the uniform measure on 
the sphere is invariant, and the generator is strongly elliptic. From the theory of irreducible Markov 
processes with compact state space, we know that the process is ergodic and thus a{z) converges for 
large z to the uniform distribution over the sphere of radius Ro- This can be used to compute the 
limit distribution of the mode powers (|aj P);,=i,...,Ar for large z, which is the uniform distribution 
over the set 

N 

Un = {{P,b=i,...,iv, P, >0,J2 Pj = ^0} ■ (4.27) 

i=i 

We carry out a more detailed analysis that is valid for any z in the next section. 

4.4.2. Independence of the change of coordinates that flatten the boundaries. The coefficients 
(|4.14p , (|4.16p and (|4.17p of the generator £ have simple expressions and are determined only by the covariance 
functions of the fluctuations v{z) and n{z) and the boundary values of the derivatives of the eigenfunctions 
(j)j{ui,^) in the unperturbed waveguide. The dispersion coefficient Kj has a more complicated expression 
(|4.2ip - (|4.23p . which involves integrals of products of the eigenfunctions and their derivatives with powers of 
^ or X — ^. These factors in ^ are present in our change of coordinates 

t{z, i) = B{z) -f [T{z) - i?(z)]l =i + e[(X^ i)^i{z) + iv{z)] , (4.28) 

so it is natural to ask if the generator C depends on the change of coordinates. We show here that this is 
not the case. 

Let F^{z,£) G ([0,00) X [0,X]) be a general change of coordinates satisfying 

for each e > 0, and converging uniformly to the identity mapping as e 0, 

sup sup |F^(z,0 -CI = 0(e), sup sup \dzF^z,£)\ ^ 0{e). (4.30) 

2>0 5G[0,X] z>OjG[0,X] 
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Note that (|4.30p is not restrictive in our context since {^^{z), v{z)) and their derivatives are uniformly bounded. 
Define the wavefield 



^, z) ^ p{uj, F%z, z) , 



(4.31) 



and decompose it into the waveguide modes, as we did for u{u!, ^, z) = p{uj, £^{z, (,), z) . We have the following 
result proved in appendix |B] 

Theorem 4.3. The amplitudes of the propagating modes of the wave field (|4.3ip converge in distribution 
as e to the same limit diffusion as in Theorem\4-.2\ 



4.4.3. The loss of coherence of the wave field. From Theorem 14.21 and the expression (|4.24p of the 

generator we get by direct calculation the following result for the mean mode amplitudes. 

Proposition 4.4. As e — > 0, E[a|(cL),z)] converges to the expectation of the limit diffusion aj{uj,z), 
given by 



E[aj (w, z)] = Oj^oi^) exp ■ 



n 



z + i 



Kj{uj) 



(4.32) 



As we remarked before, v'f^ — Tf/ is negative, so the mean mode amplitudes decay exponentially with the 



range z. Furthermore, we see from (|4T4l) and (|4T6l) that - F^''^ is the sum of terms proportional to 



.(0) 



J] 
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{d^(f)j{X)) / Pj and (9^(/>j(0)) / f5j. These terms increase with j, and they can be very large when j ^ N . 
Thus, the mean amplitudes of the high order modes decay faster in z than the ones of the low order modes. 
We return to this point in section [51 where we estimate the net attenuation of the wave field in the high 
frequency regime ^ 1. 

That the mean field decays exponentially with range implies that the wave field loses its coherence, and 
energy is transferred to its incoherent part, the fluctuations. The incoherent part of the amplitude of the 



j— th mode is — E[a^], and its intensity is given by the variance E[|a* 
if its mean amplitude is dominated by the fluctuations, that is if 



The mode is incoherent 



?e|2l 



1/2 



We know that the right hand side converges to (|4.32p as e ^ 0. We calculate next the limit of the mean 



powers 



4.4.4. Coupled power equations and equipartition of energy. As we remarked in section |4. 4. 1[ 
the mode powers |aj(i^, z)^, for j = 1, . . . , N , converge in distribution as e — )■ to the diffusion Markov 
process {Pj{uj, z))j=i,...,Ar supported in the set (|4.27p . and with infinitesimal generator Cp. We use this result 
to calculate the limit of the mean mode powers 



(cj,z) = E[Pj{uj,z)] = limE[|a^(w,z)|^ 



Proposition 4.5. As e ^ 0, E[|a^(cj, z)p] converge to Pj^\llj,z), the solution of the coupled linear 
system 



dP. 



(1) N 



dz 



0. 



(4.33) 



with initial condition P^^\uj,z = 0) = |aj^o(w)p, for j ~ 1, . . . ,N. 



Matrix F^'^' [cu) is symmetric, with rows summing to zero, by definition. Thus, we can can rewrite (|4.33p in 
vector-matrix form 



^^^=r(^)(..)p(i)(z), z>0, and p(i)(0)==Pii 



(4.34) 
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with P(i)(z) = (^P[^\ . . . , Pi^^Y and Pj^^ the vector with components |%o(a;)p, for j = 1,...,A^. The 
solution is given by the matrix exponential 

(z) = cxp [r(") {io)z\ Pj;^^ . (4.35) 

We know from (|4.14p that the off-diagonal entries in r''^-' are not negative. If we assume that they are 
strictly positive, which is equivalent to asking that the power spectral densities of ly and fi do not vanish at 
the arguments (3j — Pi, for all Z = 1, . . . , iV, we can apply the Perron- Frobenius theorem to conclude that 
zero is a simple eigenvalue of r'^-'(cj), and that all the other eigenvalues are negative, 

^N{u:)i^) < ■ ■■ < A2(CJ) <0. 

This shows that as the range z grows, the vector P(^^(z) tends to the null space of r'*^', the span of the 
vector (1, . . . , 1)^. That is to say, the mode powers converge to the uniform distribution in the set (|4.27p at 
exponential rate 

32,' 



sup 

j=l....M{uj) 



pW(^,„.)_M^ <Ce-|A^(")l^ (4.36) 
■' A* (a;) 



As z — c», we have equipartition of energy among the propagating modes. 

4.4.5. Fluctuations of the mode powers. To estimate the fluctuations of the mode powers, we use 
again Theorem 14.21 to compute the fourth order moments of the mode amplitudes: 

4')(c.,z) = lmiE[|a^^(c.,2)naf(w,z)|2] ^ E[P,{uj, z)Pi{u, z)] . 

Using the generator Cp, we get the following coupled system of ordinary differential equations for limit 
moments 



jp(2) N 

dz ^ 

n = 1 



dz 

with initial conditions 



N 

2r?P?^ + Er|:^(^n^-^^0+ErgK^-^^0' ^ > (4.37) 



" 7^ j 

2) N N 



P^?{Q)^\a,,o?M- (4.38) 



The solution of this system can be written again in terms of the exponential of the evolution matrix 

for i z£ I wp rnnrlnHp that tViit 



It is straightforward to check that the function P^^ = 1 + 5ji is a stationary solution of (|4.37p . Using 

(c) 

the positivity of TV for j ^ /, we conclude that this stationary solution is asymptotically stable, meaning 



(2) 

that the solution Pj^ ' (z) converges as z -> oo to 



1 



(2). . I N{N+l) 



2 



N{N +1) 



where = X]j=i l%,oP- This implies that the correlation of Pj{z) and Pi{z) converges to —1/{N — 1) if 
j I and to {N — 1)/{N + 1) ii j = I a.s z ^ oo. We see from the j ^ I result that if, in addition, the 
number of modes N becomes large, then the mode powers become uncorrelated. The j = I result shows 
that, whatever the number of modes N , the mode powers Pj are not statistically stable quantities in the 
limit z ^ oo, since 

Var(P,(tj,z)) ,^oo A^- 1 
E[P,(a;,z)]2 ^ iV + 1 ■ 
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5. Estimation of net diffusion. To illustrate the random boundary cumulative scattering effect over 
long ranges, we quantify in this section the diffusion coefficients F^j^ and rj^^ in the generator £ of the limit 
process. In particular, we calculate the mode-dependent net attenuation rate 



that determines the coherent (mean) amplitudes as shown in (|4.32p . The attenuation rate gives the range 
scale over which the j— th mode becomes essentially incoherent, because equations (|4.32p and (|4.35p give 

^M^iil^ «i if,»^-i. 



E 



\aj{uj,z)f 



|E[S,(w,z)]|' 



The reciprocal of the attenuation rate can therefore be interpreted as a scattering mean free path. The 
scattering mean free path is classically defined as the propagation distance beyond which the wave loses its 
coherence [2D]. Here it is mode-dependent. 

Note that the attenuation rate /Cj(w) is the sum of two terms. The first one involves the phase diffusion 
coefficient F^*^^ in the generator Cg , and determines the range scale over which the cumulative random phase 
of the amplitude aj becomes significant, thus giving exponential damping of the expected field E[aj]. The 
second term is the mode- dependent energy exchange rate 

W = (5.2) 

given by the power diffusion coefficients in the generator Cp. Each waveguide mode can be associated with 
a direction of incidence at the unperturbed boundary, and energy is exchanged between modes when they 
scatter, because of the fluctuation of the angles of incidence at the random boundaries. We can interpret 
the reciprocal of the energy exchange rate as a transport mean free path, which is classically defined as the 
distance beyond which the wave forgets its initial direction [20] . 

The third important length scale is the equipartition distance l/|A2(a;)|, defined in terms of the sec- 
ond largest eigenvalue of the matrix r'-'^-'(a;). It is the distance over which the energy becomes uniformly 
distributed over the modes, independently of the initial excitation at the source, as shown in equation (|4.36p . 

5.1. Estimates for a v^faveguide with constant wave speed. To give sharp estimates of Kj and 
J^j for J = 1, . . . , A^, we assume in this section a waveguide with constant wave speed c(^) = Co and a high 
frequency regime ^ 1 . Note from (|4.13|) that the magnitude of F^^' depends on the rate of decay of the 
power spectral densities Tl^{(3) and 'R.^{(5) with respect to the argument /3. We already made the assumption 
(|4.9|) on the decay of the power spectral densities, in order to justify the forward scattering approximation. 
In particular, we assumed that Ti-u{P) — 'R-niP) — for all (3 > 2/3jv Thus, for a given mode index j, we 
expect large terms in the sum in (|4.13p for indices I satisfying 

1/3, ~/3i\<2pN^ yV2^, (5.3) 

where we used the definition 

= -Tpy/{N + aY ~ f, j = l,...,N, and = iV + a, for a e (0, 1) . (5.4) 

A TT 

Still, it is difficult to get a precise estimate of F^-^'* given by (|4.13p . unless we make further assumptions on 
TZi, and 7?.^. For the calculations in this section we take the Gaussian covariance functions 



7^^(^) = cxp ( j and 7e^(z) = exp ( ) , (5.5) 
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and we take for convenience equal correlation lengths ii, = ifj^ ~ £ . The power spectral densities are 

ii^iP) = a^iP) = V2^^exp (^~^) , (5.6) 

and they are negligible for /3 > 3/^. Since N = [fcX/yrJ, wc sec that (j5.3p becomes 

3 27r 3 
\P] - Pi\<-< —y/2aN or equivalently, kt > ^=Vn 1 . (5.7) 
<- 2\/2a 



Thus, assumption (|4.9p amounts to having correlation lengths that are larger than the wavelength. The at- 
tenuation and exchange energy rates (|5.ip and (|5.2Pa re estimated in detailed in Appendix [C] We summarize 
the results in the following proposition, in the casqj 



iV < H < iV. (5.8) 

Proposition 5.1. The attenuation rate ICjicu) increases monotonically with the mode index j. The 
energy exchange rate J'jiui) increases monotonically with the mode index j up to the high modes of order N 
where it can decay if M 3> For the low order modes we have 

Jj{lo)X « K.j{uj)X - {M)-^/'^, j - 1 . (5.9) 
For the intermediate modes we have 

J,iuj)X ^ IC,iu;)X ^ (j7^)^ l<Cj<iV. (5.10) 

v/i-(jW 

For the high order modes we have 

J^{u)X^—, IC,{oj)X ^ klN\ j^N, (5.11) 



for kl ~ ViV, hut when ki > VN, 

Jj{uj)X < 1Cj{lo)X - HTV^ , j . (5.12) 



The results summarized in Proposition 15.11 show that scattering from the random boundaries has a 
much stronger effect on the high order modes than the low order ones. This is intuitive, because the modes 
with large index bounce more often from the boundaries. The damping rate ICj is very large, of order 
N'^ki for j ^ N, which means that the amplitudes of these modes become incoherent quickly, over scalecil 
ranges z ^ XN~'^{kt)^^ <C X. The modes with index j ^ 1 keep their coherence over ranges z = 0{X), 
because their mean amplitudes are essentially undamped K,jX <ti 1 for j ^ 1. However, the modes lose their 
coherence eventually, because the damping becomes visible at longer ranges z > X{k£y/'^. 

Note that the scattering mean free paths and the transport mean free paths are approximately the same 
for the low and intermediate index modes, but not for the high ones. The energy exchange rate for the high 
order modes may be much smaller than the attenuation rate in high frequency regimes with ki ^ \/N. These 
modes reach the boundary many times over a correlation length, at almost the same angle of incidence, so 
the exchange of energy is not efficient and it occurs only between neighboring modes. There is however a 
significant cumulative random phase in Oj for j ^ N , given by the addition of the correlated phases gathered 
over the multiple scattering events. This significant phase causes the loss of coherence of the amplitudes of 
the high order modes, the strong damping of E[aj]. 



^The case k£> N is also discussed in Appendix [Cl 
^Recall from section [4?l] that the range is actually z/e'^. 
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Note also that a direct calculatioifl of the second largest eigenvalue of r^'^\uj) gives that 

|A2H|«|rlfH|^(M)-V2. 

Thus, the equipartition distance is similar to the scattering mean free path of the first mode. This mode 
can travel longer distances than the others before it loses its coherence, but once that happens, the waves 
have entered the equipartition regime, where the energy is uniformly distributed among all the modes. The 
waves forget the initial condition at the source. 

5.2. Comparison with waveguides with internal random inhomogeneities. When we compare 
the results in Proposition 15.11 with those in [21 Chapter 20] for random waveguides with interior inhomo- 
geneities but straight boundaries, we see that even though the random amplitudes of the propagating modes 
converge to a Markov diffusion process with the same form of the generator as (|4.24l) . the net effects on 
coherence and energy exchange arc different in terms of their dependence with respect to the modes. 

Let us look in detail at the attenuation rate that determines the range scale over which the amplitudes 
of the propagating modes lose coherence. To distinguish it from (jS.ip . we denote the attenuation rate by ICj 
and the energy exchange rate by J'j, and recall from [SJ Section 20.3.1] that they are given by 



N 



J, 



3 ' 



T^3i Wj - A) 



(5.13) 



Here TZji{z) is the Fourier transform (power spectral density) of the covariance function TZji(z) of the 
stationary random processes 



X 



dx {x)(j)i {x)v{x, z) , 



the projection on the eigenfunctions of the random fluctuations v{x, z) of the wave speed. 

For our comparison we assume isotropic, stationary fluctuations with mean zero and Gaussian covariance 
function 



7^(x,z) = E [i/(x, z)i^(0, 0)] 



so the power spectral densities are 



7^,^(/3) 



Thus, (|5.13p becomes 



TtI"^ (ki)^ ( Xli \2 
p 2 \,rN ) 



\ N N J 



\ N ^ N J 



(5.14) 



J, 



' 8X (l + a/Nf ~{j/Ny 



8X 



E 

; = 1 

1^3 



{l/Nf 



\ N N } 



\ N } 



and their estimates can be obtained using the same techniques as in Appendix [C] We give here the results 
when ki satisfies (15.81). For the low order modes we have 



7r(fc£)2 



[{ktf + Nki\ ^ Nki> j 



-—^ NU> N^/^, 

8 ki 



*By direct calculation wc mean numerical calculation of the eigenvalue. We find that for TV > 20 and for kl > \/~N, 



|A2(aj)| Ri Irj^'^^cj)! with a relative ( 



that is less than 1%. 
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and for the high order modes we have 



IC-jX 



TrN{ke) 



2 r 



2s/2M 



Thus, we sec that in waveguides with internal random inhomogeneities the low order modes lose coherence 
much faster than in waveguides with random boundaries. Explicitly, coherence is lost over scaled ranges 

z <XN-^I'^ < X. 

The high order modes, with index j ~ N , lose coherence over the range scale 

z <XN~^/'^ < X. 

Moreover, the main mechanism for the loss of coherence is the exchange of energy between neighboring 
modes. That is to say, the transport mean free path is equivalent to the scattering mean free path for all the 
modes in random waveguides with interior inhomogeneities. Finally, direct (numerical) calculation shows 
that 

O {{kl)-') < ^ < O ((M)-3/2) , 
\Ji\ 

so the equipartition distance is larger by a factor of at least O (iV^/^) than the scattering or transport mean 
free path. 

6. Mixed boundary conditions. Up to now we have described in detail the wave field in waveguides 
with random boundaries and Dirichlet boundary conditions (|1.4p . In this section we extend the results to 
the case of mixed boundary conditions (|1.5p . with Dirichlet condition at x = B{z) and Neumann condition 
at a; = T{z). All permutations of Dirichlet/Ncumann conditions are of course possible, and the results can 
be readily extended. 

Similar to what we stated in section[2l the operator +uj'^c~'^{x) acting on functions in {0,X), with 
Dirichlet boundary condition at a; = and Neumann bomidary condition x ~ X, is self-adjoint in L^{0,X). 
Its spectrum is an infinite number of discrete eigenvalues Aj(a;), for j = 1,2,..., and we sort them in 
decreasing order. There is a finite number iV(w) of positive eigenvalues and an infinite number of negative 
eigenvalues. We assume as in section [2] that N{uj) ^ N is constant over the frequency band, and that 
the eigenvalues are simple. The modal wavenumbers are as before, l3j{uj) = y^\Xj(uj) \ . The eigenfunctions 
(f>j{Lu,x) are real and form an orthonormal set. 

For example, in the case of a constant wave speed c{x) = Co, we have 



A, = k' 



U - l/2)7r ' 
X 



0.(x) = ,/|sin((l^^), ,^1,2,..., (6.1) 



and the number of propagating modes is given by iV = + ij . 

6.1. Change of Coordinates. We proceed as before and straighten the boundaries using a change 
of coordinates that is slightly more complicated than before, due to the Neumann condition at x = T[z), 
where the normal is along the vector (1, — T'(z)). We let 



where 



p(i, x, z) = u(t, X{x, z), Z{x, z)) , (6.2) 



Xix,z) = X J; ^^'] (6.3) 
T(z) - B{z) 

Z{x,z) = z + xT'{z) + Q{z), Q{z)^~ [ dsT{s)T"{s). (6.4) 

Jo 
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In the new frame we get that ^ = X{x, z) G [0,X], with Dhichlet condition at ^ = 

M(t,e = 0,0-0. 

For the Neumann condition at ^ = X we use the chain rule, and rewrite 

dMt, X = T(z), z) = [d., - T'{z)d,]p{t, X = T{z), z) = 0, 

as 

d^u{t, e - ^, C = ZiT{z), z)) [~d^X + T'{z)d,X] (x ^ T{z), z) + 
d^u{t, e - ^, C = Z{T{,z), z)) [ - d.,Z + T'{z)d,Z] {x = T{z), z) ^ . 
This is the standard Neumann condition 

because 

[-d,Z + T'{z)d,Z] {x = T{z),z) = -T'(z) + T'(z) [l + T{z)T"{z) + Q'(z)] = , 

and 



(6.5) 



(6.6) 



[-d^X + T'{z)d,X] {x = T{z), z) 



X+[T'{z)Y 



^0. 



T{z)^B{z) 

Now, the method of solution is as before. Using that e is small, we obtain a perturbed wave equation 
for u, which we expand as 

Cou + £Ciu + £^C2u = 0{e'), (6.7) 

with leading order operator 

L^ = dl + dl+uj^/c\i), 

and perturbation 

Ci ^ -2{v - /i)5| + 2{X - - fi')dci - 2X{X - ^y'd^ - X{X - ^y'dc, - (6.8) 

[Xm" + ii")\d, + io\d,c-'{0) [Xfi M)e] . 

6.2. Coupled Amplitude Equations. We proceed as in section [X^ We find that the complex mode 
amplitudes satisfy p.l6p - p.l7p with C, instead of z, where the ^-dependent coupling coefficients are 



with 



C^C)^eCf{0+e^cf{C)+O{e'), 



1 



1 



2 / daX~ 04>jdi^i 



1 



X 



X 



2vW 



X 



di{x-0<l>Ai 



(6.9) 
(6.10) 

(6.11) 
(6.12) 
(6.13) 
(6.14) 



_____ l'2\ 

and coefficients c^j-; and df^^.ji defined by (|3.48l) and p.50p . Similar formulas hold for Cj; (C). 

In the following we neglect for simplicity the evanescent modes, which only add a dispersive (frequency 
dependent phase modulation) net effect in the problem. These modes can be included in the analysis using 
a similar method to that in section [ 
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6.3. The Coupled Mode Diffusion Process. As we have done in section SI we study under the 
forward scattering approximation tlie long range limit of the forward propagating mode amplitudes. 

First, we give a lemma which shows that the description of the wave field in the variables (x, z) or Q 
is asymptotically equivalent. 

Lemma 6.1. We have uniformly in x 

X ^x, — ^ — X 0, Z ^x, — ^ 2 ^ E[z/'(0)^]z in probability. 

Proof. The convergence oi X to x is evident from definitions (|6.3|) and p.2p . Moreover, ()6.4|) gives 

and integrating by parts and using the assumption that the fluctuations vanish at z = 0, we get 

^ J) - J = =^ [(- - (J) - (J) (J)] + ■ 

The first term of the right-hand side is of order e and the second term converges almost surely to E[;^'(0)^]2; 
which gives the result. □ 

The diffusion limit is similar to that in section WM and the result is as follows. 

Proposition 6.2. The complex mode amplitudes {a^j{ui,C,))j=i^,,,^N converge in distribution as £ — >■ 
to a diffusion Markov process process (aj(aj, C))j=i^....Ar. Writing 

the infinitesimal generator of the limiting diffusion process 

C = Cp+Ce 

is of the form (^T77p, but with different expressions of the coefficients given below. 

(c) 

The coefficients Tj^' in Cp are given by 

r;.f(w)=^^(/3,-/30Q^,,i+^M(ft-A)OLi if^'^'' (6-15) 

where 

= c„j; + - /?,) - (A ~ [eu,ji + f.MPi ~ ft)] 

X ' ^ 



c{xY 



Ml 



ct>j{X)UX), (6.16) 



The coefficients in Cg are similar, 

rf^\u;) = ^^(o)g^,,, + Uf^ml^i Vj, z , (6.17) 



and 



rlf (^) = 7.,,7g^ .7 + l^..,iQln (6.18) 



J' 

with 7j^j7 and 7^^; defined by (j4.18p . 

We find again that these effective coupling coefficients depend only on the behaviors of the mode profiles 
close to the boundaries. In the case of Dirichlet boundary conditions, the mode coupling coefficient rj;''(cj) 
depends on the value of d^(j)jd^4>i at the boundaries. In the case of Neumann boundary conditions, the mode 
coupling coefficient r^;''(aj) depends on the value of (f)j{X)4>i{X). 

Given the generator, the analysis of the loss of coherence, and of the mode powers is the same as in 
sections 14.4. 3M. 4. 51 
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7. Summary. In this paper we obtain a rigorous quantitative analysis of wave propagation in two 
dimensional waveguides with random and stationary fluctuations of the boundaries, and either Dirichlet or 
Neumann boundary conditions. The fluctuations arc small, of order e, but their effect becomes significant 
over long ranges zje^ . We carry the analysis in three main steps: First, we change coordinates to straighten 
the boundaries and obtain a wave equation with random coefficients. Second, we decompose the wave 
field in propagating and evanescent modes, with random complex amplitudes satisfying a random system of 
coupled differential equations. We analyze the evanescent modes and show how to obtain a closed system 
of differential equations for the amplitudes of the propagating modes. In the third step we analyze the 
amplitudes of the propagating modes in the long range limit, and show that the result is independent of the 
particular choice of the change of the coordinates in the first step. The limit process is a Markov diffusion 
with coefficients in the infinitesimal generator given explicitly in terms of the covariance of the boundary 
fiuctuations. Using this limit process, we quantify mode by mode the loss of coherence and the exchange 
(diffusion) of energy between modes induced by scattering at the random boundaries. 

The long range diffusion limit is similar to that in random waveguides with interior inhomogeneities 
and straight boundaries, in the sense that the infinitesimal generators have the same form. However, the 
net scattering effects are very different. Wc quantify them explicitly in a high frequency regime, in the 
case of a constant wave speed, and compare the results with those in waveguides with interior random 
inhomogeneities. In particular, we estimate three important length scales: the scattering mean free path, 
the transport mean free path and the equipartition distance. The first two give the distances over which 
the waves lose their coherence and forget their direction, respectively. The last is the distance over which 
the cumulative scattering distributes the energy uniformly among the modes, independently of the initial 
conditions at the source. 

We obtain that in waveguides with random boundaries the lower order modes have a longer scattering 
mean free path, which is comparable to the transport mean free path and, remarkably to the equipartition 
distance. The high order modes lose coherence rapidly, they have a short scattering mean free path, and 
do not exchange energy efficiently with the other modes. They also have a transport mean free path that 
exceeds the scattering mean free path. In contrast, in waveguides with interior random inhomogeneities, all 
the modes lose their coherence over much shorter distances than in waveguides with random boundaries. 
Moreover, the main mechanism of loss of coherence is the exchange of energy with the nearby modes, so 
the scattering mean free paths and the transport mean free paths are similar for all the modes. Finally, the 
equipartition distance is much longer than the distance over which all the modes lose their coherence. 

These results are useful in applications such as imaging with remote sensor arrays. Understanding 
how the waves lose coherence is essential in imaging, because it allows the design of robust methodologies 
that produce reliable, statistically stable images in noisy environments that we model mathematically with 
random processes. An example of a statistically stable imaging approach guided by the theory in random 
waveguides with internal inhomogeneities is in [3]. 
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Appendix A. Proof of Lemma 13. 11 The proof given here relics on explicit estimates of the series in 
(|3.33p . obtained under the assumption that the background speed is constant c(^) = Co- Wc rewrite (|3.33p 
as 



(A.l) 



with linear integral operators ^'i and ^'2 defined component wise by 



Tb \ iMji~d,Q'^i){z + s)diiw,z + s)e 




(A.2) 




(A.3) 
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The coefficients have the exphcit form 



= h Hz) - Kz)] (f )' + ""^'\^"^"H 6,i + (1 - 5,i) K(z) - ^"(z)^ 



(1 - [1 - (-1)'+^] + 0{e), (A.4) 

Q5,(z) = [v\z) - ^,'{z)] 5,1 + (1 - 5,i) W{z) - ^,'{z)] -i^ - 

(1 - 5,iy{z)-^ [1 - (-1)'+^] + 0(e). (A.5) 

Let £^(Z;L^(M)) be the space of square summable sequences of L^(R) functions with hnear weights, 
equipped with the norm 



'illL2(K))^ 



1/2 



We prove that ^ : t{{Z\ L^{R)) t{{Z\ L^{R)) is bounded. The proof consists of three steps: 

Step 1: Let T be an auxihary operator acting on sequences v = \vi\i^i^ defined component wise by 

[T.], = E = E (tTI + 7^) 4 ((-^ ^-'^ * 1 + -0 * y)^. + - V,). 

This operator is essentially the sum of two discrete Hilbert transforms, satisfying the sharp estimates 

||f * y||£2 < 7r||t;||f2. 

Therefore, the operator T is bounded as 

\\Tv\\,. < {1/2 + n) EII^jII^?- (A-6) 



Step 2: Let v{z) = be a sequence of functions in K and define the operator 

Q : ^ il(Z;L'{R)), [Qv],iz) = [Tt;], * e-ftl«l(z) l{,>jv}, (A.7) 

where 



Using Young's inequality 

||[Q«],|U.(K) = ||[rt;], *e-ft-l^l|U.(K) < ||[rt;],|U.(M)||e-ft-l^l|Ui(R) - A ||[Ti;],|U.(r), (A.9) 
we obtain from (|A.6|) - (|A.9p that HQH < (1 + 2tt)/C{uj), because 

EOii[Q^].iuW'^cf^ Eii[r^biiW) = 7^ / Eit^^'^b-wi'rf^ 

<7vJ;w(V2 + 7r)^/ Eb-^.-(^)l^'^^- '^cLV^' EOll-.-|l^-(«))'- 
^ -^^ jGZ ^ ^ jm 
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This estimate applies to the operator ^'2. Indeed, let us express ^'2 in terms of the operator Q using 
(D and (1X51) . 

['^2vUz) = - ^^')v,) * e-ftl^l(z)l{,>^.} - 2[Qi,' vi^z) + 2{-iy[Qiy'{-iy m],{z). (A.ll) 

That the sum in 4*2 is for / > A'^ is easily fixed by using the truncation vi = vi l{;>jv}- Thus, using estimate 
(jA.lOp for the last two terms, we obtain 

l!*2t^||£f < -^(^ (llMlki>~(R) + lklki>~(R)) \\v\\ii 



step 3: It remains to show that the operator is bounded. We see from (jA.2[) . (jA.4[) and (|A.5[) that 
for any j > N 



2,; 2 



1 



Pj 



where ^£'2 is just like the operator ^1/2, with the driving process {u' , fi') replaced by its derivative {i>" , fi" 
Using again Young's inequality, we have 

mivUL^m ll(''"'^)^jHi'w + -cRll[*2«]jllL-(R)- 

Now multiply by j and use the triangle inequality to obtain that ^'i is bounded, 

2tt^ , (5 + 87r) 



(lkllM/2,=c + 11^11,4,^2,00) 



Appendix B. Independence of the change of coordinates. We begin the proof of Theorem 
with the observation that 

where is the inverse of £^ , meaning that -w and u are related by composition of the change of coordinate 
mappings. Clearly, the composition inherits the uniform convergence property 



sup sup =0(£)- 

2>o ?e[o,x] 



(B.l) 



For the sake of simplicity we neglect the evanescent modes in the proof, but they can be added using 
the techniques described in section [231 Using the propagating mode representation of ^, z), 



N 



N 



= ^(/)i(a;,Owi(a;,z) +^4>i{ui, ^, z)ui{uj, z), 



(B.2) 



1=1 



1=1 



where we let 







{i^^-\z, F%z, 0) - 9^0, {u;, s i'^-\z, F%z, 0) + (1 - s) ds. 



But we can also carry out the mode decomposition directly on w and obtain 

N 



w{uj,£,,z) = ^(f>i{uj,^)wi{uj,z), 



(B.3) 



1=1 
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because the number of propagating modes N and the eigenfunctions in the ideal waveguide are independent 
of the change of coordinates. Here wi{uj,z) are the amphtudes of the propagating modes of w. Equating 
identities (|B.2p and (|B.3p . multiplying by and integrating in [0,X] we conclude that 



Wj{uj,z) = Uj{uj,z) + ^Qj(a;,z)uz(w,z), (B.4) 



where we introduced the random processes, 

ci,{lo,z)^ r <j>,{uj,0 [' d^cl,i{uj,st'-\z,F'{z,0) + {^~s)i) {z, {z, 0) ~ dsd^. 

Jo Jo 

In addition, differentiating equation (jB.4p in 2, we have 

N 

dzWj (w, z) = dzUj (w, z) + ^ dzCij (w, z)ui {oj,z) + cij (w, z)dzUi {uj,z). (B.5) 

Now, let us recall from the definition of the forward and backward propagating modes that 

il3jUj{uj, z) + dzUj{oj, z) = 2i^/]Tj aj{uj, z)e''^'"'. 
We conclude from (|B.4j) and (|B.5p that 

+ ^-^''^"^ fc(c.,.)e-(^^-ft)- + &K^,^)e-'^^-+^')-) , (B.6) 

i=i V Pi 

where {a"'(aj, z)}j=i jv arc the amplitudes of the forward propagating modes of w{io^^,z). A similar 

equation holds for the backward propagating mode amplitudes {bj{uj, z)}j=i n- 

The processes cij{uj, z) can be bounded as (|4.30p 

,i^'?'3^»r'f™Pl^y('^'^)l} - ™P sup |a5(/)z(a;,0|} X 

i<J,l<N z>0 i<3,l<N ^e[o,X] Ce[o,A'] 

sup sup \t'-\z,F%z,0)-^\^O{e). (B.7) 
2>o ?e[o,x] 

For the processes dzCij{u;, z) we find a similar estimate. Indeed, note that 

a, [a^^, (c, s i^^-\z, F%z, 0) + (1 - s) (^''"'(^: ^'(^, e)) - e)] = 

-A, ^'■-i(z,F^(z,0) + {l~s)C)s d,[t^-\z,F^z,m {e''-\z,F'{z,0) - + 

+ (1 ~ a.[r'-i(z,^^^(z,c))]. 

A direct calculation shows that 

X{F^{z,0-e^l{z)) 



A(l + £v{z)) - eyL{z)_ 
(d^F^z, - e^i'iz))iXil + eu{z)) - e^(.T)) - {F^{z, ^ ~ e^i{z)) e {u'{z) - ^i' {z)) 



{X{l + ev{z))-e^i{x)f 
Hence, using condition (|4.30p for dzF^{z,£) 

sup sup \dz[t^~\z,F'{z,0)\\<C{\\v\\w^,^,M\w^,^)e. 

^>0 ?6[0,X] 
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Therefore, 

max {sup|5zQj(cj,z)|} < X max {A; sup \(j>j{uj,^)\ sup |0;(a;, ^)|} O(e^) + 

X max { sup sup la^f/)^^., 01} 0(e). (B.8) 

l<j.l<N ^g[o,x] «6[0,X] 

Let a™(Lu,z) and be the vectors containing the forward and backward propagating mode am- 

phtudes and define the joint process of propagating mode amphtudes X^{z) = {a^{u!, z),b'^{uj, z))'^ . Let 
us the long range scaled process be X^'''"{z) — X™(z/e^). Equation (|B.6p implies that 



X-j^iz) = X-Jz) + M, [u;,C[u;,-j,d.C[L,,-),-) X^z), (B.9) 

where C(aj, z) :~ (q^ (oj, z))j./=i^....jv and dzC(j-L>,z) := (92Qj(w, z))j.(=i....,Ar. The subscript e in the matrix 
Me(-) denotes the fact that this matrix depends explicitly on e and. due to estimates (jB.7|) and (jB.8|) . we 
have 

sup ||Me(c^, C{uj, z), a,C(w, z), z)||oo = 0(e). (B.IO) 

Let us prove then, that the processes X^''^{z) and X^{z) converge in distribution to the same diffusion 
limit. Denote by Q{Xo, L) the 2A''-dimensional cube with center Xq and side L. The probability that 
Xfj'"'{z) is in this cube can be calculated using (|B.9|) . 



K--(z) € Q(Xo, L)] = / dP- (a;, 4) 



{xe(I+M,(C,9,C,z))-iQ(a!!o,L)} 



z 



a;,C,a,C,^ . (B.ll) 



Here P*" (a;, z) is the probability distribution of the process X™ (z) and P (a;, C, 9zC, z) is the joint probability 
distribution of the processes {Xi^{z), C(a;, z), 9zC(w, z)). We can take the inverse of I + Me(C, S^C, z) by 
()B.10p . The same estimate (|B.10|) also implies that for every 5 > Q there exists eo such that for e < eo, 

{x e Q{xo, (1 - 5)L)} C {a; e (I + M,(C, 9,C, z))-iQ(a;o, L)} C {a; G Q(a;o, (1 + 5)L)}. (B.12) 

Denote the diffusion limits by 

X^{z) = lim X;(z), X»(z) = Imi X^-«'(z). 
We conclude from (|B.11[) and (jB.12[) that for any 5 > 0, 

P[X^(z) e g(Xo, (1 - <5)L)] < PK'(^) e g(Xo, i)] < P[X„(z) e g(Xo, (l + 5)L)]. 
Sending (5 — )■ 0, we have that for any arbitrary cube g(a;o, L) 

P[X<,(z) e g(Xo,L)] =P[X-(z) e g(Xo,L)]. 

This proves that the limit processes have the same distribution and therefore, the same generator. 

Appendix C. Proof of Proposition [5711 Recall the expression (|2.3p of the wavenumbers. The first 
term in (|5.ip follows from (|4.16p : 



- (J)' [«.(0) + «,.(0)] jf^^^, » ^llFT^F^ ^ 



p(0) 



It increases monotonically with j , with minimum value 



(0) (27r)3/2 M , , 

rl?-^^^«i, (C.2) 
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and maximum value 

(27r)3/2 



2aX 

The second term in (jS.ip . which is in (|5.2p . follows from (14.131) . (j5.6p and (j5.4 



If < j/N < 1, then we can estimate (|C.4I) by using the fact that the main contribution to the sum in I 
comes from the terms with indices I close to j, provided that k£ is larger than N^^^ and smaller than N. 
We find after the change of index I ^ j + q: 

^■'^ ^ X((7V + a)2-j2)iV^^ 
Interpreting this sum as the Riemann sum of a continuous integral, we get 



- r'^^'luj) ~ — — —— / e ds = . ^ ' ' — (C.5) 

^^^^ X((7V + a)2-j2)y_^ XV(iV + «)2-j-2 ^ ^ 

By comparing with (jC.l|) we find that the coefficient — T^^'' (w) is larger than rj°' when satisfies \fN <^ 

k£<t:N. 

To be complete, note that: 

- If H - N, then -T^f^iuj) is larger than rj°' if and only if j/N < (1 + (fc£/A^)2)- V2. 

- If fc€ is larger than iV, then the main contribution to the sum in I comes only from one or two terms with 
indices I — j ±1, and it becomes exponentially small in (fc£)2/A^2_ j.^ ^j^ese conditions — r^^''(a;) becomes 
smaller than r^^-* . 

For J ^ 1 we can estimate (|C.4p again by interpreting the sum over I as a Riemann sum approximation 
of an integral that we can estimate using the Laplace perturbation method. Explicitly, for = 1 we have 

A Jo V 1 - S2 



We approximate the integral with Watson's lemma [2] Section 6.4], after changing variables ^ = (1 — \/l — ^2)2 
and obtaining that 

Watson's lemma gives 

1 



and therefore by (|C.6[) and ([F 



^-^(i-^Y ^ r(3/4)2V4 

X/T^ ^ (fc£)3/2 
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By comparing with (|C.2[) we find that the coefficient —T^ii{uj) is larger than T^ii ■ 

For j ^ N only the terms with I ~ N contribute to the sum in ()C.4[) . If k£ ^ VN , then we find that 

^^^^ ^V^T^ 2C{a)MX' 

up to a constant C(a) that depends only on a. By comparing with (jC.Sp we can see that it is of the same 
order as fI^L. If kl ^ -y/iV, then we find that 



which is very small because the exponential term is exponentially small in {klY /N . In these conditions 
— F^^^(a;) is smaller than F^j^. 
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